dd <- c("C:/xuncao/mcmc_net_polity/2002try 2") 
setwd(dd)

dat<-dget("data_wd_nomissing_polity2002")

Xs<-na.omit(dat[[3]])
Xs[, 4]<-log(Xs[, 4])
Xs[, 5]<-log(Xs[, 5])

Xr<-Xs

Yd<-dat[[1]][rownames(Xs), rownames(Xs)]
Yd[is.na(Yd)]<-0 # be consistent with Xd
diag(Yd)<-NA
Yd<-log(Yd+1)

Xd<-dat[[2]][rownames(Xs), rownames(Xs), ] 


####
####

OUT=read.table("C:/xuncao/mcmc_net_polity/2001try 2/OUT",header=T) 
U=read.table("C:/xuncao/mcmc_net_polity/2001try 2/U")
V=read.table("C:/xuncao/mcmc_net_polity/2001try 2/V")
 
nss=dim(OUT)[1]
n=dim(U)[1]/nss
k<-3
nburn<-100
aa<-colnames(Yd)

PU<-array(dim=c(n,k,nss))
for(i in 1:nss) { PU[,,i]<-as.matrix(U[ ((i-1)*n+1):(i*n) ,])  }
PU<-PU[,,-(1:nburn)];rownames(PU)=aa;     #drop burn in
PV<-array(dim=c(n,k,nss))
for(i in 1:nss) { PV[,,i]<-as.matrix(V[ ((i-1)*n+1):(i*n) ,])  }

# drop burn in ##
PV<-PV[,,-(1:nburn)]; rownames(PV)=aa; 
OUT=OUT[-(1:nburn),]

 
#find posterior mean of Z %*% t(Z)
UTV<-matrix(0,n,n)
for(i in 1:dim(PU)[3] ) { ## Do not change k here
    UTV<-UTV+PU[,,i]%*%t(PV[,,i]) 
}
UTV<-UTV/dim(PU)[3]  ## Do not change to k here
tmp<-svd(UTV)
 
u2=tmp$u[,1:k];rownames(u2)=aa;
v2=tmp$v[,1:k];rownames(v2)=aa;
# u2[aa1,]=u1;v2[aa1,]=v1;d2=tmp$d[1:k];



#Xd=array(dim=c(length(Xs),length(Xs),2));rownames(Xd)=aa
#rownames(Xs)=rownames(Y)
#Xr<-Xs

#    Xd[,,1]=u2%*%diag(d2)%*%t(v2)
#    Xd[,,2]=xd




####

# Xd<-Xd[,,-2] # this get rid of common language ...
Xd[,,2]<-UTV # plug in the UTV from the latent space of the previous year.  



source("gbme.asym.r")
n<-dim(Yd)[1]
Y<-Yd

# call the main estimation package, gbme
gbme(Y=Yd,Xd=Xd,Xs=Xs,Xr=Xr,k=3,NS=100000,awrite=T,bwrite=T, oround=8) 
